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Abstract 

Improved aerodynamic mathematical models, for 
use in aircraft simulation or flight control design, are 
required when representing nonlinear unsteady 
aerodynamics. A key limitation of conventional 
aerodynamic models is the inability to map frequency 
and amplitude dependent data into the equations of 
motion directly. In an effort to obtain a more general 
formulation of the aerodynamic model, researchers 
have been led to a parallel requirement for more 
general testing methods. Testing for a more 
comprehensive model can lead to a very time 
consuming number of tests especially if traditional 
single frequency harmonic testing is attempted. This 
paper presents an alternative to traditional single- 
frequency forced-oscillation testing by utilizing 
Schroeder sweeps to efficiently obtain the frequency 
response of the unsteady aerodynamic model. 
Schroeder inputs provide signals with a flat power 
spectrum over a specified frequency band. For 
comparison, experimental results using the traditional 
single-frequency inputs are also considered. A method 
for data analysis to determine an adequate unsteady 
aerodynamic model is presented. Discussion of 
associated issues that arise during this type of analysis 
and comparison of results using traditional single 
frequency analysis are provided. 
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Nomenclature 

A, B, C numerator transfer function coefficients 
a, b, indicial function parameters 

c mean aerodynamic chord, m 

C L lift coefficient 

Cm pitching moment coefficient 

J cost function 

k non-dimensional frequency, k= 10 / /V 

/ characteristic length, / = c / 2 

N number of frequencies 

PSD power spectral density 

q pitch rate, rad/sec 

s Laplace transform variable 

t time, sec 

V airspeed, m/sec 

v measurement noise 

Z output measurement vector 

a angle of attack, rad 

0 unknown parameter vector 

a 2 variance 

t dummy integration variable 

to angular frequency, rad/sec 

superscripts 

estimate 

subscripts 

L lift force 

m pitching moment 
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Introduction 

Conventional aerodynamic math models for use in 
aircraft simulation or flight control design have become 
increasingly deficient as aircraft maneuvering 
capabilities have extended the flight envelope beyond 
conventional boundaries. In these new flight regimes, 
especially cases where aircraft maneuver at high 
angular rate or high angle of attack (a), nonlinear 
unsteady aerodynamic effects associated with rigid 
body responses have been well documented. The need 
to improve the aerodynamic mathematical model in 
these cases is generally acknowledged 1 3 . Although the 
aerodynamic modeling problem is readily identified 
with military fighter aircraft, it is also an issue for 
civilian transport aircraft. Airliner loss-of-control 
accidents due to adverse weather or system failures 
place these aircraft in nonlinear flight regimes that are 
not well modeled. 

The deficiency in most aircraft simulations and 
models used in control design is the use of a traditional 
formulation for the aerodynamic model 6 . This 
formulation assumes that the aerodynamic forces and 
moments can be represented by a differentiable 
function and therefore expanded in a Taylor Series 
expansion with only first order terms (stability and 
control derivatives) retained. Although this 
mathematical structure is adequate for benign portions 
of the flight envelope, it is completely incapable of 
modeling aircraft responses when nonlinear unsteady 
effects occur. Information about the aerodynamics can 
be obtained from static and dynamic wind tunnel 
testing. 

A key limitation of the conventional model 
structure is inability to map frequency and amplitude 
dependent dynamic data into the equations of motion 
directly. In addition, this formulation only 
accommodates the measurement of combined stability 
derivatives during wind tunnel oscillatory tests 
preventing separate estimates of damping and unsteady 
terms. Because the conventional combined terms are 
frequency, amplitude, and angle of attack dependent, 
simulation engineers are forced to find ad hoc methods 
to incorporate unsteady effects into the equations of 
motion 5 . 

In an effort to obtain a more general formulation 
of the aerodynamic model, researchers have been led to 
a parallel requirement for more general testing 
methods. Tobak 7 and others have suggested testing 
methods that allow measurement of the basic 
parameters of a general aerodynamic model. Testing 


for a more comprehensive model can lead to a very 
large number of tests, especially if traditional single 
frequency harmonic testing is attempted. This paper 
offers an alternative to traditional single-frequency 
testing by utilizing Schroeder 8 sweeps to obtain the 
unsteady aerodynamic model. Schroeder inputs 
provide signals with a flat power spectrum over a 
specified frequency band with low peak-to-peak input 
amplitudes. Amplitudes are controlled by properly 
phasing each harmonic. 

The modeling of aircraft aerodynamics in planar 
one degree -of-freedom motion in pitch is identical to 
that in Refs. [12] and [13]. For data analysis the least 
squares and maximum likelihood estimates were used. 
The methodology is tested on wind tunnel data from a 
10% scale model of an F-16XL configuration. The 
results are compared with those obtained from 
traditional static and single frequency data in Ref. [12] 
and large amplitude oscillatory data in Ref. [13]. 

Modeling and Identification 

Results from wind tunnel forced-oscillation tests 
for given amplitude and mean angle of attack show that 
the resulting combination of stability derivatives 
depends on the frequency of the oscillations. This 
dependency contradicts the basic assumption that 
stability derivatives are time invariant. The effect of 
frequency on the aerodynamic parameters was 
explained by Goman 7, 10 et al„ at TsAGI, and by 
Klein 11, 12 et a!., at NASA LaRC, by formulating the 
linear aerodynamic model equations with unsteady 
terms. 

As an example, lift is considered in the form 

t 

C L (0 = C L (0) + J C La 0 - T)a(T)dT 

0 ( 1 ) 
/ r 

+—) C L q (r-TMr) dT 

0 

where C L (I ) and C L ( t. ) are the indicial functions, 

C L (0) is the initial value of C L , I is the characteristic 
length, and V is the airspeed. Two assumptions were 
adopted to allow simplification of the model used in the 
analysis of measured data: a) the effect of q(t) on the 
lift can be neglected and b) the indicial function 
C L ( t. ) can be expressed as 
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C La (t) = C La (oo)- ae - b S. (2) 



H 


The simplified model, which takes into accomit 
increments witli respect to steady state conditions, can 
be expressed as 

C L (0 = C La (oo)a(t) + -Lc Lq (oo )q(t) 

r (3) 

— aj e~ h[ i,T> a(T)dT 
o 

where C L (°°) and C, (°°) are the rates of change 

with a and q evaluated in steady flow. 

Applying the Laplace transform to Eq. (3), the 
expression for the lift coefficient is obtained as 


B=C L Joo)- a + b 1 ^C Lq (o o) (8) 

C=b y C La { oo). 

Because Eq. (5) is nonlinear in the parameters, the 
estimation represents a nonlinear estimation problem. 
The initial values of parameters for this technique were 
obtained from a linear regression using the cost 
function 

N 2 

J{6) = C L (m + i<Oj) + ( Aa> 2 -C-iBcOj)a(j ) | (9) 
y'=i 

Experiment 




As 2 +Bs + C 
s+b 1 


a(s) 


(4) 


where s is the Laplace transform parameter. 


The problem addressed in this paper is the 
identification of the model given in (4). The parameters 
in (4) can be estimated by applying the maximum 
likelihood principle in frequency domain to the 
dynamic wind tunnel data. The model used in the 
estimation procedure has the form 


C L (co) 


- Am 2 +C + iBa) 
+ ico 


a(co) 


(5) 


-U) = C L (j) + v(j), j = 1,2,..., N (6) 


where C L (w) and a( w) are the Fourier transforms of 
C L (it) and a(t) , v( j) is the measurement noise 
assumed to be a Gaussian random complex sequence 
with zero mean and variance a’, N is the number of 
frequencies at which the transformed input output data 
are known, and co is the angular frequency. 

The maximum likelihood estimator minimizes the 
negative logarithm of the likelihood function 

0 = min \-k\L{Z N -,9,<7 2 * )\ (7) 

e,a 2 1 J 

where Z N = [z(l), z(2), . . ., z(N)] is a vector of output 
measurements and 9 = [A, B, C, bi] is the vector of 
unknown parameters. The parameters in (3) were 
obtained by solving the equations 


A drawing of the 10% scale F-16XL model is 
shown in Figure 1 . Dynamic tests were conducted in 
the NASA Langley 12-Foot Low-Speed Wind Tunnel. 
For these tests the model was mounted on a dynamic 
test rig through a six -component strain-gauge balance. 
The dynamic test rig is a computer controlled 
hydraulically actuated system, which was sting- 
mounted on a C-strut support system. The mounting 
arrangement rotated the model about the reference 
center of gravity location of 0.558 c . The tests were 
conducted at a dynamic pressure of 192 Pa (4 psf) 
resulting in the flow velocity of 17.52 m/sec (58 fps) 
and a Reynolds number of 10 6 based on the mean 
aerodynamic chord. 

Oscillatory data were created using Schroeder 
sweeps in a as an input at 1 3 mean values of angle of 
attack, Oo, and an amplitude, a A = 5 degrees. The range 
of a-mean values was from «,, = 0 to 65 degrees. Data 
were sampled at 100 Hz with an analog filter at 10 Hz. 
Tests were repeated ten times at each angle of attack 
and then an average signal was formed using the 
ensemble data. The ensemble-averaged data was used 
for data analysis. 

Example time histories of angle of attack, lift and 
pitching moment coefficients, a, C L , and Cm, at oto = 

36 degrees is given in Figure 2. These plots show 
displacements relative to starting values at a®. The 
harmonic content of the angle of attack is shown in 
Figure 3 as a function of reduced frequency, k. Figure 3 
indicates a flat spectrum for a frequency range of 0.03 
Hz to 2.0 Hz which corresponds to a range of non- 
dimensional frequencies k = 0.004 to 0.268. For the 
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analysis, time histories of a and the longitudinal 
aerodynamic coefficients were transformed to the 
frequency domain using a Discrete Fourier Transform 
(DFT) algorithm 14 . This algorithm allowed the 
transform to be performed over the frequency range of 
0 to 2 Hz. 

In addition, the traditional in-phase and out-of 
phase components of aerodynamic coefficients were 
computed from transformed data. Expressing 
transformed angle of attack and lift coefficients as 
a(co) = u x + iv 2 and C L (co) = u 2 + iv 2 , it can be shown 
that the transfer function 


C L (co) 

a(a>) 


U 2 +iv 2 
u t + i\\ 


= U +iV 


( 10 ) 


Then the in-phase and out-phase components are, 
respectively, 


C 


= C^-k 2 C L2l =U and 


C L =C L +C Ld =-V. 

A A La k 


Both Fourier components are obtained from measured 
data and plotted against angle of attack for three values 
of reduced frequency. These plots are presented in 
Figure 4. 

Results and Discussion 

The estimation method used in this study is 
different from that used in References [11-13] because 
it provides estimates of all unknown parameters at 
different values of a. The previous methods in 
References [11] and [12] assumed that the tune 
constant, 1 /b , , did not change with a. In Reference [13] 
it was assumed that initial estimates of C L (°°) were 

available from small-amplitude oscillation data to 
initialize estimation of the unsteady term in Eq. (3). 

The parameter estimation method outlined in this 
paper was applied to longitudinal (planar) low- 
amplitude forced-oscillation data. The method was only 
applied to cases in the range of angle of attack between 
a = 20 degrees and a = 65 degrees. This range of angle 
of attack is where unsteady effects primarily occur for 
this model. The estimated parameters in Eq. (3) and 
their 2-sigma bounds are shown in Figures 5 and 6 for 
the lift and pitching moment equations, respectively. 
The figures also show alternate estimates of the 
parameters when possible. 


Figure 5 presents estimates of unknown 
parameters in model Eq. (3). The estimation process 
failed for a < 25 degrees because of the lack of 
unsteady information in the measurement data. This 
identifiability issue is corroborated in Figure 4, where 
the variation of the in-phase and out-of-phase 
components with frequency is greatly reduced. Below a 
= 25 degrees there is virtually no variation with 
frequency. 

Estimates of C L (°o) show very good agreement 

with results obtained from static measurements. Static 
data estimates were obtained by performing a simple 
finite difference calculation on C L ( a). Estimates of 
C L (°°) , on the other hand, did not show good 

agreement with values in Ref. [12] and error bounds 
were relatively large. This result is not too surprising 
since C L (°°) is usually a very small contributor to lift 

and therefore difficult to estimate. Estimates of 
parameter a, the gain on the unsteady term, agreed vety 
well with two other independent techniques in 
References [12] and [13]. Again the error bounds are 
large where the information content of the measured 
data is substantially reduced. Estimates of bj generally 
follow the same variation with a as that found in 
Reference [13]. The match with Reference [12] is only 
approximate since the parameter was assumed to be 
constant in that study to minimize identifiability issues. 

In Figure 6, for a < 45 degrees, estimates 
of C "',, (°°) show very good agreement with estimates 

obtained from simple finite difference calculations on 
the static data. Increasing error bounds for estimates in 
the high and low a ranges occur as before due to low 
information content. Poor signal-to-noise ratio in the 
pitching moment channel also contributed to poor 
identifiability, particularly in the high a range. In 
Figure 2 an apparently more noisy response is seen in 
the pitching moment channel at a = 36 degrees. Signal- 
to-noise ratio deteriorated further with increasing angle 
of attack. This problem was aggravated by the limited 
capability of the dynamic rig to follow low-amplitude 
oscillation commands. The rig was originally intended 
for large amplitude and relatively fast oscillations, hi 
this low-amplitude experiment the Schroeder sweep 
required a number of small oscillations as small as one 
degree, particularly at the beginning of the sweep (see 
Figure 2). The pitch channel was particularly sensitive 
to this problem and for the a = 46 degrees case 
convergence for the estimation algorithm was poor. 
Further investigation is required to fully explain the 
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pitching moment response. Excluding the a = 63 
degrees case, estimates of C m (°°) match very well 

with values given in Ref. [13] and the 2-sigma bounds 
were relatively small. Estimates of parameter a, for the 
pitching moment equation, also agreed very well with 
the independent technique in Ref. [13] for a < 45 
degrees. For a > 45 degrees, estimates differed from 
Ref. [13] and large error bounds were obtained that 
reflect the poor information content and low signal-to- 
noise ratio. In addition, parameter a tended to be very 
small in the low and high a ranges making the 
estimation more difficult. Estimates of b show a trend 
toward zero or small values for a above 45 degrees, 
thus indicating a trend toward neutral stability. With 
parameter a reduced the unsteady effect in C m is limited 
for a greater than 45 degrees. As mentioned previously, 
the estimates of bi from Ref. [12] assumed no 
dependence on a. Estimates of the same parameters in 
Ref. [13] didn’t reveal any dependencies on a. Below 
45 degrees, estimates of bi were in the same range as 
estimates from both Refs. [12] and [13]. 

Figures 7 and 8 show the frequency response of 
the lift and pitching moment coefficients, respectively, 
during forced oscillations at a = 36 degrees. These 
figures show the relatively large measurement noise 
(frequency domain) for this study. Amplitude ratio and 
phase plots for the lift coefficient are relatively flat at 
the high frequency range (2 Hz). The pitching moment 
amplitude ratio is small and slightly increasing with 
frequency. 

An identifiability issue of concern for the transfer 
function model structure given in Eq. (10) is pole/zero 
cancellation. Figure 9 shows the pole and two zeros as 
a function of a for the lift transfer function. The 2 nd 
zero shown as a diamond in the upper plot is relatively 
large and has a wide frequency separation from the 
other transfer function roots in the lower angle of attack 
range. The pole and 1 st zero are fairly close in 
magnitude and become almost equal at very low angle 
of attack. The lower plot highlights the pole and 1 st 
zero to show how the two approach cancellation at a = 
20 degrees. At that point the values are within 5% of 
each other creating a more difficult estimation problem. 
As the pole and 1 st zero values approach equality the 
model structure in Eq. (10) becomes inadequate and the 
estimation algorithm fails. 

Concluding Remarks 

This study is part of an ongoing effort at NASA 
Langley to develop a more general formulation of the 


aerodynamic model for aircraft that includes nonlinear 
unsteady aerodynamics. The general structure proposed 
for the aerodynamic model has the form 

^ ^ static L dynamic ^ unsteady 

L unsteady = ~ (0 

where 

C* L = —b(a)C* L -l- a(a)a 

Similar equations can be written for drag and pitching 
moment equations. This structure allows easy 
interpretation of the model parameters by retaining 
conventional stability and control derivatives for static 
and dynamic terms. Unsteady terms are obtained by 
solving a first order differential equation with a- 
dependent coefficients. This approach offers a 
straightforward model for estimation and simulation. 

Modeling unsteady aerodynamics demands 
substantial testing in wind tunnels and in turn a parallel 
requirement for more general and efficient test 
methods. In support of these goals this study 
demonstrated an alternative to single-frequency 
harmonic testing using Schroeder sweeps. A transfer 
function model that includes unsteady aerodynamics 
was developed. A frequency domain method for data 
analysis that can be used to estimate the unsteady 
aircraft models was presented. Identifiability issues 
associated with this methodology and with this 
particular experiment were discussed. 

The methodology in this study compared well 
with previous studies that used different techniques to 
obtain the aerodynamic model for the F-16XL. Both 
present and past techniques separated the static, rotary, 
and unsteady terms and both modeled the unsteady 
term as an indicial function. However, the previous 
techniques made different assumptions about a priori 
information. The previous methods assumed that some 
of the parameters were known. The present approach 
has the advantage that all unknown parameters in the 
model are estimated at once for various angle of attack. 

The previous techniques only used conventional 
single-harmonic, constant-amplitude inputs requiring a 
substantially larger number of tests. This study 
demonstrated the use of Schroeder sweeps that 
provided wide-band. flat-spectra inputs with low peak- 
to-peak input amplitudes. The maximum amplitude for 
the Schroeder sweep is constrained by the dynamic 
rig’s capabilities and the bandwidth required for the 
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experiment. The Schroeder algorithm prescribes 
smaller amplitude cycles within the sweep that improve 
the input from an experiment design point of view, 
especially for nonlinear or input-amplitude dependent 
systems. In general, the Schroeder sweep is an excellent 
frequency-rich input for forced oscillation testing 
however the dynamic rig must be able to handle the 
variety of input amplitudes and frequencies. When test 
rig capabilities are limited, an alternate choice may be a 
traditional frequency sweep with limited amplitude 
variation. 

Identifiability issues for the unsteady model, 
discussed in this paper, were associated with either a 
lack of information content or an inadequate model 
structure. Degradation of information content in system 
response was found in three cases. First, low levels of 
frequency dependence were found for in-phase and out- 
of-phase coefficients in the low and high a ranges. 
Second, low signal-to-noise ratios were found in 
pitching moment responses, especially in the high a 
range or low amplitude cases. Third, when parameter 
values were very small causing little contribution to the 
overall response. Inadequate model structure was 
highlighted in this study primarily for tire problem of 
pole-zero cancellation at very low angles of attack. 
Other contributors to a mismatch between model 
structure and system response are suspected, i.e., the 
dynamic rig, designed for large amplitude motion, may 
have introduced responses not associated with the 
aerodynamics. This was especially a problem in the 
pitching moment channel at mid to high angles of 
attack and for small amplitude motion. In this study 
responses in the pitch channel did not completely 
follow the expected harmonic responses. Only minor 
distortion was expected since even at 5 degrees 
amplitude some nonlinear behavior can occur in the 
pitch channel. 
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Figure 1. Three-view drawing of 10% F-16XL model. 
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Figure 2. Tune histories of angle of attack, lift, and 
pitch-moment coefficients for 0(0=36 degrees. 
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Figure 3. Harmonic content of transformed a input for 
oLq = 36 degrees. 
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Figure 4. Variation of in-phase and out-of-phase 
components of lift coefficient with angle of atack for 
three different reduced frequencies. 
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Figure 5. Estimated parameters and their 2a 
confidence bounds for lift coefficient. 


Figure 6. Estimated parameters and then 2a 
confidence bounds for pitch-moment coefficient. 
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Figure 7. Comparison of measured and computed 
frequency response, Q/<x 
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Figure 9. Pole and zero variation with angle of attack 
for lift coefficient transfer function. 




Figure 8. Comparison of measured and computed 
frequency response, C,Ja. 
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